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ABSTRACT 

During their main sequence lifetimes, the majority of all Galactic Disc field stars must 
endure at least one stellar intruder passing within a few hundred AU. Mounting obser- 
vations of planet-star separations near or beyond this distance suggest that these close 
encounters may fundamentally shape currently-observed orbital architectures and hence 
obscure primordial orbital features. We consider the commonly-occurring fast close en- 
counters of two single-planet systems in the Galactic Disc, and investigate the resulting 
change in the planetary eccentricity and semimajor axis. We derive explicit 4-body an- 
alytical limits for these variations and present numerical cross-sections which can be 
applied to localized regions of the Galaxy. We find that each wide-orbit planet has a few 
percent chance of escape and an eccentricity that will typically change by at least 0.1 
due to these encounters. The orbital properties established at formation of millions of 
tight-orbit Milky Way exoplanets are likely to be disrupted. 

Key words: planets and satellites: dynamical evolution and stability - planet-star 
interactions - stars: kinematics and dynamics - Galaxy: kinematics and dynamics - 
Galaxy: structure - celestial mechanics 



1 INTRODUCTION 

After leaving their birth clusters, most stars undertake a po- 
tentially harrowing multi-Gyr journey through the Galactic 
Disc. The stars are continuously perturbed by global Galactic 
phenomena and are periodically nudged by individual stel- 
lar encounters. Occasionally, an encounter is close enough 
to cause major disruption to any planets orbiting in the ap- 
proaching systems. The currently observed exoplanet popu- 
lation may be shaped in part by these encounters. 



1.1 Typical Closest Encounter Distances 

Using simple a r gume nts (e.g. from Pg. 3 of 
iBinnev fc Tremaind 120081 ). one can crudely estimate an 
upper bound for the typical encounter distance, r enc , over 
a main sequence lifetime. If n denotes the space density of 
stars in the Galactic Disc, and u ran is the random velocity 
of stars, then r enc w (47TOii ra ntMs) _1 ' /2 , where tyts is the 
main sequence lifetime. This estimate is conservatively 
large because gravitational focusing is not included. We 
can estimate tMs through simulations from the SSE stellar 
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evolution code (|Hurlev et al.l [20001. Doing so yields Fig. [TJ 
which plots the closest encounter main sequence distance 
as a function of progenitor mass from IMq — 2Mq, which 
represents a common range of exoplanet host masses. The 
majority of stars drawn from a standard stellar initial mass 
function (see e.g. IParravano et alj l201lf ) will have masses 
under 1 Mq, further suggesting that the typical encounter 
separations in Fig. [1] represent overestimates. The solid and 
dashed lines represent Solar and very low (l/200th of Solar) 
metallicities, respectively. The metallicity of a star helps 
dictate its main sequence lifetime, and hence the expected 
close encounter distance. The plot partially illustrates that 
differences in the metallicity of stars have little (indirect) 
effect on the close encounter distance. 



The figure demonstrates that the majority of all stars 
will suffer a close encounter of just a few hundred AU for a 
reasonable range of n and i» ra n values. Even in sparse envi- 
ronments, like the Solar neighborhood (with n ~ 0.1 pc -3 ), 
Sun-like stars will approach one another at least once within 
a few hundred AU. This e stimate corroborates the roug h esti- 
mate of 500 AU given by Zakamsk a fc Tremaind (|2004T ). who 
consider only a 5 Gyr encounter timescale. 
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Figure 1. Upper bound estimates for typical main-sequence clos- 
est encounter distances, r enc , between exosystems in the field of 
the Galactic Disc. Solar metallicity stars and very low metallicity 
stars are represented by solid and dashed lines, respectively. The 
random velocity of stars is « ran and the space density of stars is n. 
Fiducial values of n (0.6pc — 3 ) and i> ra n (50km/s) are represented 
by the black curves. Because the majority of Galactic Disc stars 
are less massive than IMq, they will have r enc values less than 
hundreds of AU. 



1.2 Stellar Encounter Orientations with Respect to 
Galactic Centre 

Given that close encounters within hundreds of AU will typ- 
ically occur, we can now attempt to characterize the orienta- 
tions of the collisions with respec t to the Galactic Centre. As 
outlined by iQuillen et al.l (|201ll 'l. the distribution of veloci- 
ties in the Galactic Disc is affected by a multitude of factors. 
Pote ntial perturbers in c lude Galactic Lind blad resonances 
(e.g. I Yuan fc Kuolll997l ; iLepine et al]|201ll '), stellar streams 
from past mergers and interactions with satellit e subhaloes 
(e.g. iBekki fc Freeman! 120031 ; iGomez et all l20ld), and tran- 
sient spiral density waves fe.g. |Pe Simone et alj|2004 ). Stel- 
lar velocities may also be highly dependent on t he phase and 
patte rn speed of the Milky Way's spiral arms (|Antoia et alj 
l201lf ). suggesting drastic differences in the velocity distribu- 
tion in different regions. These factors might help explain 
why the velocity components of the stars in th e Solar neigh- 
borhood are neither isotr opic nor Gaussian l|Binnev et alj 
l200d : lNakai ima et ai"1l2010h . Generally, the orbits of Disc stars 
are modulated vertically and epicyclically (e.g. pgs. 164-166, 
iBinnev fc Tremaine| l2008 ) , and may undergo significant ra- 
dial migration fe.g. ISchoenrichll201lf ). Further, the amplitude 
of the epicyclic and vertic al oscillations are of the sa me order 
of mag nitude (e.g., Pg. 18. lBinnev fc TremainefeoOct ), and are 
orders of magnitude longer than then the physical radii of the 
stars themselves. Therefore, we should expect that stars suf- 
fer close encounters with each other at random orientations 
with respect to the Galactic Centre. 



1.3 Planetary Orbit Orientations with Respect to 
the Galactic Disc 

Now we assess whether the planes of the planetary orbits 
should have a preferential orientation to the Galactic Disc. 



The severe misalignment of the Sola r System's invariable 
plane with the Gala ctic plane at w 60° (|Huang fc Waddl 19661 ; 
iDuncan et al.|[T987l ) foreshadows the likely answer. Observa- 
tions constrain the distribution of exoplanetary orbital planes 
poorly because most extrasolar planets have been discovered 
by the Doppler radial velocity technique, which alone does 
not provide any information about line-of-sight inclinations. 
Similarly, the stellar rotational axis orientation - which is 
suggestive of planetary orbit orientation - of the vast major- 
ity of non-exoplanet host stars is unknown. Howe v er, in cases 
where this info r matio n has been obtained, lAbtl (|200ll ) and 
iHowe fc Clarkd ( 2009) find that the these axes are orientated 
randomly. For exoplanet-host stars that harbour transiting 
planets, we do have line-of-sight inclination information. Ac- 
cording to the Exoplanet Data ExploreiQ, as of 15 January, 
2012, there are 141 transiting exoplanets. At most, the orbital 
plane of any of these planets is misaligned with our line-of- 
sight by ~ 13.4°. However, the median misalignment angle is 
just w 2.79° and the standard deviation is ~ 2.73°. 

Therefore, effectively we observe transiting planets edge- 
on, and the locations of these planets on the sky might suggest 
a relation between the Galactic plane and planetary orbital 
planes. In Fig. [2] we plot the declination versus right ascension 
of the host stars of these 141 transiting planets from this 
database. In order to help assuage the strong observational 
bias in the plot, the plot markers are colored and shaped 
according to the planet names, which are often indicative of 
the program or collaboration who first discovered the planet. 
For example, the planets with names containing "Kepler" or 
"KOI" (Kepler Object of Interest) are all clustered in the 
same region on the plot. This is due to the the fixed field 
the Kepler space mission is observing. The plot definitively 
illustrates that observed planetary orbital planes are known 
to encompass a wide range of orientations with respect to the 
Galactic Disc. 

These considerations lead us to treat close encounters 
between two planetary systems in arbitrary directions with 
orbital planes that are arbitrarily aligned with each other. 
However, we must sensibly restrict the vast phase space of 
these encounters. We do so first by reviewing some published 
literature related to this topic. 



1.4 Extending Previous Scattering Studies 

The three-body problem which includes a star-planet pair 
experiencing a perturbation from an intruder star has been 
the subject of several studies, and is well-characterized 
in many regions of phase space. Most studies, however, 
treat the se interactions in t h e context of cluster en- 
counters dHeggie fc Rasiol Il99rj; iDavies fc Sigurdssonl boOll ; 



iFregeau et al.M2006l ; ISpurzem et all 120091) . which typically 
have a higher n, smaller u ran , and much shorter lifc- 
time (tens of Myr) t han in the field. An exception is 
IZakamska fc Tremainel (|2004r ) , who do consider perturbations 
in the field from a stellar intruder, but on a multi-planet sys- 
tem. They treat the perturbation as a superposition of three- 
body interactions, neglecting the contribution to the potential 
from the planets. Also, they treat the velocity vector of the 



1 See the Exoplanet Data Explorer at |http:/ /e xoplanets. org/ 
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Figure 2. Approximate line-of-sight exoplanetary orbital plane 
orientations. Plotted are the spatial coordinates of stars which 
host transiting planets. All data is taken from the Exoplanet Data 
Explorer, as of 15 January, 2012. Plot markers are determined 
based on whether the orbiting planet's name includes "WASP" 
(blue filled circles), "HAT-P" (downward-pointing brown filled tri- 
angles), "OGLE" (hollow gray squares) "Kepler" or "KOI" (hollow 
red diamonds), "TrES" (purple filled diamonds) or "XO" (upward- 
pointing orange filled triangles). Other transiting planets are given 
by filled black squares. The plot demonstrates that planetary or- 
bital planes are known to encompass a wide range of alignments 
with respect to the Galactic Disc. 



intruder and the orbital plane as coplanar. Among their sev- 
eral useful results are i) about 10% of all stars experience 
close encounters within 200 AU, ii) planetary eccentricities 
may be excited up to 0.1 in the field, and iii) the extent of 
the excitation is strongly dependent on system size and phase. 

Here, we provide a multi-tiered extension to that 
work. First, we consider the potential of all four bod- 
ies in the close encounter of two one-planet systems, as 
most Milky Way s tars are now thought to have planets 
(|Cassan et al.l |2012| ). Previous studies of the 4-body prob- 
lem often consider the more general case of the interaction 
of two stellar binaries (Mikkola 1984: Hut 1993: Bacon et al 



e s W 

19961 : iHeggiel I2OO0I: iGiersz fc Spurzeml 120031: iFregeau et al 



2004 : iPfahl fc Muterspaughl 120061 : ISweatmanl 120071 ) or a 
planet-less intruder perturbing a mult i-planet exosystem 
dMalmberg et al.ll201ll : iBolev et al"1l2012h . However, none of 
these studies consider the close encounter of two single-planet 
systems. 

Second, because field encounters are fast, we develop an 
analytical method based on impulses that can determine the 
change in orbital parameters without resorting to numerical 
simulations. We consider two extremes in phase for our anal- 
ysis, although the method can in principle be generalized to 
arbitrary phases, and even arbitrary numbers of planets. 

Third, we do perform numerical simulations, here specif- 
ically for the purpose of obtaining normalized cross sections. 
These quantities then enable one to determine the overall 
rate of encounters and eccentricity excitations over a main 
sequence lifetime in localized patches of the Milky Way. As 
already argued earlier, we consider encounters of all mutual 



orientations, independent of their locations with respect to 
the Galactic Centre. 



1.5 Plan for Paper 

We outline some of the key quantities in the hyperbolic 4- 
body problem in Section 2 before our analytical (Section 3 
and the Appendix) and numerical (Section 4) explorations. 
Of particular note are the two specific orientations we model 
without numerical integrations (Sections 13. 2. 31 and 13. 2. 4|l and 
the eccentricity excitation frequencies arising from our nu- 
merical integrations (Section I4.3[l . In Section 5, we interpret 
the results. Section 6 discusses related topics, and Section 7 
provides a short conclusion. 



1.5.1 Variables used 

Table 1 delineates the variables applied throughout this pa- 
per. The subscript k takes the values "1" and "2" and is 
used to describe the planet and star belonging to the differ- 
ent planetary systems taking part in the encounter. Primed 
and double-primed values are explained in the text where nec- 
essary. 



2 4-BODY PROBLEM SETUP 

Consider a planet with mass M p \ orbiting a Galactic Disc 
star with mass M s i, and an independently evolving planet 
with mass M P 2 orbiting a different Disc star with mass M„2. 
Initially, assume the distance between the systems (denoted 
"1" and "2") is infinity. Each planetary orbit is described by 
the planet's semimajor axis, at, and eccentricity, e^, where 
k = 1 or 2 depending on the planet. At t = 0, the orbital 
parameters are denoted by an additional subscript, "0" . 

As argued in Section 1, the systems may approach each 
other at any orientation, and the relative orientation of the 
planetary orbital planes is also unconstrained. Now consider 
the plane in which the stars approach each other, and fix the 
reference frame on M S 2 - System #1 will approach System #2 
such that M s i will be traveling at a velocity Voo with an im- 
pact parameter b. Because Voo > the stars will approach 
each other on approximate hyperbolic orbits, approximate 
because of the presence of the planets. Denote the reduced 
mass of the 2-body hyperbolic system as /i = G (M 3 i + M S 2). 
We treat values of Voo, b and fi, as well as all four individ- 
ual masses and ato, for k = 1, 2, as given, known quantities 
throughout the paper. Further, eto = always. 



2.1 Key Orbital Parameters 

The total energy of the system is equal to = —fx/2ah, 

where ah is the (negative) hyperbolic semimajor axis. Hence 



Oh 



JL 

v 2 

1 .X' 



(1) 



The total angular momentum of the system is equal to 6 Voo, 
which can be related to the hyperbolic eccentricity, e/, > 1, 
such that 



2 , , rvz 

Sh = H 5- 



(2) 
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Table 1. Variables Used in this Paper 



Variable Explanation 



a h 

(*) 
fcO 

„(*) 

a kf 

(*) 

X 



a 



a 



a 
b 

CJC 
CJC,1 

cjc , 2 

b< c > 

cjc, 3 



bpl p 2 
b p ls2 

b s lp2 

6 (c) 

"stat,< 

6 (c) 
°stat,> 



7 
<5 

Sk 
e (c) 

cxt,max 

e (c) 

cxt,mm 

eh 
M 

e kf 

J*) 
e x 

E h 

e 

V 

M pk 

M sk 

Mtot 

A 4 

n 

A" 

JV 

Q 

r enc 
^start 

RAND 
a 

& norm 
£cnc 

^integrate 
*MS 

T fc 
T 

fran 

Voo 

^circ, A; 

^circ, fcO 

K;irc,0 

K;rit 

^pcri 

| AVI | 



Hyperbolic semimajor axis for a star 

Initial semimajor axis for planet k in the far (* = /) and close (* = c) cases 
Final semimajor axis for planet k in the far (* = /) and close (* = c) cases 
Contribution to Planet #2's semimajor axis variation due to Planet #1 alone in the far (* 

Number of planetary orbital periods to numerically integrate before the close encounter 
Impact parameter of both stars 

far case impact parameter value at which a planet escapes 

Maximum close case impact parameter value separating planetary escape from boundedness 
Middle close case impact parameter value separating planetary escape from boundedness 
Minimum close case impact parameter value separating planetary escape from boundedness 

Maximum impact parameter used in the numerical simulations 

Impact parameter which causes a planet-planet collision 

Impact parameter of both planets 

Impact parameter of Planet #1 and Star #2 

Impact parameter of Star #1 and Planet #2 

close case lower impact parameter value at which there is no net perturbation on the planets 
close case upper impact parameter value at which there is no net perturbation on the planets 

Factor by which (aio + a2o) is multiplied to obtain q for the numerical simulations 
Fraction of the innermost planetary orbit used as a numerical integration timestep bound 
Dimensionless planet/star mass ratio for each system when both are physically equivalent 
Dimensionlcss planet/star mass ratio for system k 



■■ f) and close (* = c) cases 



close case local eccentricity maximum, for (b. 



.)<*> 



(c) (c) 

close case local eccentricity minimum, for feV 1 < b < feV 2 
Hyperbolic eccentricity of a star 

Final eccentricity for planet k in the far (* = /) and close (* = c) cases 

Contribution to Planet #2's eccentricity variation due to Planet #1 alone in the far (* = /) and close (* 

Hyperbolic anomaly of a star 
Dimensionless ratio equal to aio/a20 
Dimensionless ratio equal to Voo/V^ r it 
Mass of planet k 
Mass of star k 

Total mass of the 4-body system 

Sum of both stellar masses, times the Gravitational Constant 
Space density of stars 
Number of experiments 

Number of times over a main sequence lifetime that |Aei| > T occurs 
Pericenter of the star-star hyperbolic orbit 

Typical closest encounter distance for two stars in the Galactic Disc 
Separation used to initialize numerical integrations 
Low-discrepancy quasi-random Niederreiter number between and 1 
Cross section 
Normalized cross section 

Timescale of close encounter between both planetary systems 

Numerical integration timescale 

Main Sequence lifetime 

Orbital period of planet k about star k 

Given extent of an eccentricity perturbation 

Random stellar velocity 

Velocity of Star #1 with respect to Star #2 at an infinite separation 

Circular velocity of planet k about star k 

Circular velocity of planet k about star k assuming M pk = 

Circular velocity of either planet for equal planetary masses and semimajor axes 
Velocity at which the total energy of the 4-body system equals zero 
Pericenter velocity of the star-star hyperbolic orbit 

Magnitude of the velocity kick perpendicular to the direction of motion 



: c) cases 
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The pericenter, q > 0, of the star-star hyperbolic orbit is 
/. . b 2 V£ 



\a h \ (e h - 1) 



V- 



1 



(3) 



2.2 Velocity Comparisons 

The total energy of a 2-body system with a nonzero rela- 
tive velocity is positive. The critical velocity of the four body 
system, V cr it, for which the total system en ergy is zero and 
ionization is possible is l|Fregeau et al.ll2004T ): 



Vcrit 



' GM tot 


( M s iM pl 


M s2 M p2 \ 


(M sl + M p i) (M s2 + M p2 ) 




«20 / 



where M to t = M s i +M p i + M s2 + M P 2. We plot typical values 
of V^rit in Fig. showing that Krit is nearly 23 times lower 
for two Mj planets and two Mq stars than for four Mq stars, 
where Mj is the mass of Jupiter. Hence, comparison of typical 
stellar velocities in the Galactic Disc (w lOkm/s - lOOkm/s) 
implies that one-planet systems are moving too fast to ionize 
all four bodies through encounters regardless of the values of 
aio and 020- 

Now we can compare the circular velocity of a planet with 
respect to its parent star, V C ir C ,k, to typical values of Voo. We 
have 



V:irc.fc — 29. T9 



Msk±M I k\ /1AU\ 
M Q J V Ofco J 



(5) 



Therefore, for wide orbit planets and typical Disc veloci- 
ties, Voo 3> 14i rCj fc. However, for planets on tight orbits, 
the velocities are comparable. Further, we denote V c i rc ,fco as 
the circular velocity of planet k when M p u = (such that 

Vcirc,fc ~ V:irc,fco)- 

The fastest velocity achieved in a hyperbolic orbit is at 
the pericenter of that orbit. The pericenter velocity Vperi, is 
related to Vx> through 



Vperi — 



Ijl. 


~e h + 1" 


4 \a h \ 


_e h - 1 



+ 1 



1 + ^ 



(6) 



which is always greater than Voo and becomes infinite as b 
0. Voo represents the minimum velocity of the orbit. 



where k" indicates the planet with the smaller orbital period. 
As demonstrated by Eq. (|7|), the impulse approximation is 
well-suited for wide orbits due to the resulting low value of 
Kirc,fc". In the impulse regime, the planets do not progress 
in their orbits around their parent stars during the encounter 
(i.e., the mean anomaly is approximated as stationary). 

The impulse approximation allows us to isolate and es- 
timate analytically the planets' mutual perturbations dur- 
ing the encounter. For simplicity, let us treat both planets 
on circular orbits. By symmetry, in the impulse approxi- 
mation the only net perturbation is perpendicular to the 
velocity vector of the perturber. For ease of reference to 



IZakamska fc Tremaind (|2004T ). we also take both planetary 
(4) systems to be coplanar with each other and with the per- 
turber's velocity vector. We will be estimating the pertur- 
bations on Planet #2. By symmetry, the perturbations on 
Planet #1 will yield the same change in orbital parameters. 

3.1 General Case 

3.1.1 Total Perturbations on the Passing Star 

First, let us estimat e the perturbations o n Sta r #2 due to 
Star #1. Pg. 422 of iBinnev fc Tremaind (|l987h shows that 
the imparted velocity kick is 



2bVl _ 2 
e h . 



(8) 



Planet #1 will also kick Star #2. The effective impact 
parameter between Planet #1 and Star #2, b p i S 2, will depend 
on the planet's position during the encounter. We have, 



IAV 



J_|pls2 



2M p i6pi,2V2, 



1 + 



u P ls2 ''oo 



G (M p i + M s2 y \ G 2 (M p i + M s 



.(9) 



3.1.2 Total Perturbations on the Passing Planet 

Similarly to the impulse imparted on Star #2 by Planet #1, 
the impulse imparted by Star #1 on Planet #2 is: 

-1 



IAV, 



slp2 



26 slp2 



1 + 



h 2 v 4 



G{M sl + M p2 ) y ' G 2 (M sl +M P 2) 2 
The impulse on Planet #2 from Planet #1 is: 

26 p ip2 V^ 



IAV 



_L|plp2 — 



G(Mj 



pi 



h 2 V 4 
°plp2 v 00 



G 2 (M pl + M p2 ) 2 



.(10) 



.(11) 



3 IMPULSE ANALYTICS 

Although we must resort to numerical simulations to fully 
explore the 4-body problem consisting of two planet-star sys- 
tems, here we investigate how this cases of this problem may 
be solved analytically in the impulse regime. As suggested by 
Eq. |[SJ|, perturbations on wide orbit planets due to passing 
plan etary systems may be tr e ated in the impulse approxima- 
tion. IZakamska fc Tremaind ([2004) claim that this assump- 
tion holds for their planetless intruder if the stellar perturber 
is fast and if the planetary period is much longer than the 
characteristic timescale of the encounter, t cn c ~ b/V per i. This 
condition is analogous here to 



T, 



2tt- 



Vt,, 



b V c j rc fc/ 



» 1, 



(7) 



3.1.3 Net Perturbations on the Passing Planet 

Therefore, Planet #2 experiences a net velocity kick relative 
to its parent star of 



|AVj_|p 2 = I A VI | s i 8 2 + |AVj_| p i a2 

- (|AVx| s ip2 + |AVx|pip 2 ) (12) 

Equipped with Eqs. (|51)- (I12I) . we can insert the se velocity 
kicks into the formalism of Jackson fc Wvattl (|2012h fsee their 
Eqs. 1-2] Q in order to determine the resulting variation in a 2 



2 The assumption under which these equations are derived is that 
the impulse is instantaneous, which is equivalent to our Eq. (7). 
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Figure 3. The critical velocity as a function of semimajor axes for two IMq stars and two lMj planets (left axis) and four IMq stars (right 
axis). Hence, given typical field velocities, unlike pairs of binary stars in the Galactic Disc, almost never will two single-planet exosystems 
be completely ionized. 



and e2. Let us denote the post-encounter values of ai and C2 
as a,2 f and e 2 f (recall 620 = 0). Then 



zero. In other specific cases of interest, these assumptions may 
be lifted and the more general results rederived in a similar 
manner as below. 



O20 

a 2f 

&2f 



IAVj. 



P2 



I AVI 



lei 



p2 



Kirc,2 



(13) 
(14) 



3.2 Specific Example 

We wish to relate aif and e2/ to b, Voo and M s i, M S 2, M p i, 
and M P 2 in an analytically tractable manner. Thus, we will 
focus on two specific cases of interest, as illustrated in Fig. [4] 
In the first case, which we denote by far, both planets are 
the furthest possible distance from each other as the systems 
pass each other; both stars are in-between the planets. Here, 
the value of b may be any value from to 00. In the second 
case, which we denote by close, the direction of the vectors 
from each star to its child planet are pointing towards each 
other. Here, a value of b that we denote 6 m i n will cause both 
planets to collide. For 6 < b m i n , the orbits will overlap. Figure 
|4]shows a cartoon of the encounters at pericenter for different 
cases. 

We have derived analytical formulae for the critical 
points of the motion, asymptotic limits, and the individual 
contribution to the perturbations from the planets alone. All 
these formulae are presented and explained in the Appendix 
in order to help retain the focus of the reader here. In this 
section, we provide just the most important results. 



3.2.1 Analytic Simplification 

In order to obtain compact, understandable formulae, for the 
remainder of Section 3, we assume aio = 020, M p i = M P 2 
and M s i = M a 2 such that both systems are equivalent ex- 
cept for their labels. Define 5 = M p2 /M S 2 and 14irc,o as the 
circular velocity of either planet assuming the planet mass is 



3.2.2 Fiducial Sample 

In order to provide tangible numbers that accompany the 
analytics and resulting plots, we concurrently consider fidu- 
cial values of M s2 = M Q , M p2 = Mj, a 20 = 1000 AU and 
Voo = 30km s _1 « 1000K:rit unless otherwise indicatecjf]. 
These values give V c i IC ,2 ~ 1.3km s _1 such that the impulse 
approximation is valid as long as b <^ 1.45 x 10 AU (Eq.Q. 



3.2.3 The far Case 

Even though the planets are at opposition to each other in 
the far case, planetary ejection will occur when the stars 
have a close-enough encounter, when b ^ , . Alternatively, 

for 6 > bJl, the orbital parameter evolution is 



p (/) 
-2f 



AS) 



2 (2a 20 + b) 



020 + b 



('211 



P 



„(/) 2 

-2/ 



(15) 



(16) 



Note that 



Also, 



,(/) 



and a^t 



a 2 o as 6 — > 00, as expected. 



a 2 y cannot decrease due to the close encounter. 



Equations (|15[) - (|16l ) show that as long as the planetary 
mass is neligible compared to the stellar mass, the planetary 
contribution is also negligible everywhere in the far case pa- 
rameter space. Nevertheless, we quantify this contribution in 
the Appendix. Note that e^J — > and a^ — > a 2 o as b 



3 Our choice of fiducial semimajor axis helps us demonstrate all of 
the regimes of interest for realistic close enc ounter distances (Fig . 
[Tj| and known exoplanet separations (e.g . I Goldman" et al.l l2010l : 
iKuzuhara et al.ll201ll ; iLuhman et~alll201li) . 
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FAR case 




CLOSE case fori >b 




CLOSE case for 6 <b 



a 20 




Figure 4. Cartoon of different close approach configurations mod- 
elled by impulses. The larger objects are stars and the smaller ob- 
jects are planets. Different colors denote the two different systems. 



as expected. Also, a^y cannot decrease due to the close en- 
counter. 

Figure [5] illustrates these properties. Depending on Voo, 

Planet 4t2 will be ejected when b is within a few tens or 
( f) 

hundreds of AU; fo^jc i s marked on the upper axis of the left 
panel for the slowest Voo- Planets which remain bound after 
surviving a passing star at b as 500 AU expand their orbits 
by tens to hundreds of AU and stretch their orbits through 
eccentricity increases of at least 0.1. 



FAR Case Impulse Approximation 




1000 2000 3000 4000 5000 
b/AU 



FAR Case Impulse Approximation 




1000 2000 3000 4000 5000 
b/AU 

Figure 5. Eccentricity and semimajor axis variation of Planet 
#2 in the far case. The variation monotonically increases as b 
decreases from infinity to when the planet is ejected. The 

eccentricity and semimajor axis always increase due to interactions 
in the Far case. 



3.2.4 The close Case 

Now let us consider the opposite limit, where the position vec- 
tors from each star to their orbiting planet point towards each 
other. The resulting orbital parameter evolution is a more 
complicated function of b. 

In particular, there are two separate ranges of b in which a 
planet will escape: i) < b < 6^ 3 , when both stars are in be- 
tween both planets and the stars are close to each other, and 
ii) ^eje 2 < & < ^eie i ' wnen one star nearly collides with one 
planet. Further, there is one region, (6^at,<) < b < ( & stit,>) , 
where the planet-planet interaction becomes important. Ad- 
ditionally, there are two local extrema: i) in between the es- 
cape regions, the perturbations are minimized at e^ t min , and 
ii) for b well beyond b m i u , the perturbations are maximized at 
e (c) 

cxt ,max ' 

All of the physical features mentioned above and illus- 
trated in both Figs. [6] and IA2I (which can be used as guides 
for the location of the critical points) can be reproduced with 
a compact analytical form for e^t as a function of impact 
parameter. We remind the reader that this quantity, among 
several others, are derived in the Appendix: 
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Figure 6. Eccentricity and semimajor axis variation of Planet #2 in the close case. At b = AU, the stars collide. At b m 1000 AU, Planet 
#2 collides with Star #1. At b ?s 2000 AU, Planet #2 collides with Planet #1. The extreme points present in the panels are explained 
in Eqs. < l A15b -< TA25t ■ The explicit functional form of e!, c ' is given in Eqs. 1171 - Both panels demonstrate that planets will experience 
major disruption and potentially ejection if their orbits cross. 



e (c) f 
such that 



(?) 



Kir 



%G 2 M 2 8l + Z 4 + Z 5 



Z(iZ7 



(17) 



Z 6 = 2GM sl + V£ (a 20 - b) 
Z 7 = 2GMsi + Vl (2a 20 - b) 



(20) 
(21) 



Z4, = 2GM sl V^(6a 20 -b(3-28)) (18) 
Z 5 = (4a 2 20 + b 2 (1 - 25) - 2ba 20 (2 - S)) (19) 



where a^f is derived from e£j in the usual way (Eq. I16[) 
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3.3 Consequences 



4.1 Numerical Simulation Setup 



The analytics show that a planet's eccentricity can be raised 
to any value due to a realistic close encounter. Even in the 
limiting case where both planets are furthest from each other 
during the encounter, if the stars endure a close enough ap- 
proach, then the planets will be ejected. In the other extreme, 
measurable eccentricity excitation can occur over a wide, re- 
alistic range of impact parameters. When a star crosses in 
between another star and planet, the eccentricity excitation 
is at least a tenth, but is likely many tenths. When any two 
bodies narrowly miss each other, planetary escape may occur. 
However, there are locations at which this net perturbation 
is zero; this range of locations increases along with planetary 
mass. 

This analytic exploration helps us to gauge expectations 
for the outcomes of numerical simulations, and perhaps more 
importantly, provides an explanation for some of the trends 
seen in the outputs of our numerical simulations. We describe 
these simulations in the following section. 



4 NUMERICAL CROSS SECTIONS 

We now calculate cross sections of various encounter out- 
comes via numerical scattering experiments. Cross sections 
of this type, introdu ced in stellar dynamical research by 
iHut fe Bahcail (|l983h . represent the effective surface area for 
some outcome of a scattering event between involving stel- 
lar or planetary systems. Coupled with a velocity distribu- 
tion and density of systems, the cross section yields a total 
outcome frequency. By suitably setting up random stellar en- 
counters and performing many experiments, the probabilistic 
outcome of these potentially chaotic encounters can be ob- 
tained. 

Here, we are interested in the frequency of planetary sys- 
tems whose planets have eccentricities that are perturbed by 
a particular amount, T. The cross section is a function of T, 
e = aio/020, and 77 = Voo/Krit such that one example is: 



a(\ Aei| > T, e, 77) = nb 



2 

max 



JV(|Aei| >T,e,t ? ) 

A'total 



(22) 



where JVtotai represents the total number of experiments and 
AT(|Aei| > T, e, 77) represents the number of experiments with 
a given e and 77 that yield |Aei| > T. One may compute er- 
rors in cr(| Aei| > T, e, 77) by using Gaussian counting statis- 
tics (|Hut fc Bahcall l983). Doing so gives error bars which are 
equal to the RHS of Eq. (32} divided by ■ s /NQAe^\ > T,e,7]). 
In order to create a scale-free cross section - for wider appli- 
cations - a can be normalized as: 



cr norm (|Aei| > T,e, 77) = 



q(|Aeij > T, e, 77) 

tv (a 10 + a 2 o) 
tr(|Aei|>T,e,»j) 



1 + e 



(23) 



Our goal is to compute values of both cr norm (| Aeij > T, e, 77) 
and cr n0 rm(| Ae2| > T,e, 77) for different values of T, e and 77. 
Doing so requires a careful numerical setup. 



First, we must set up initial conditions such that the initial 
separation of the stars is finite, and then select this finite 
separation. We also must chose a sufficiently representative 
range of b small enough to not be computationally prohibitive 
but large enough to encompass all the regimes in, for example, 
Fig. HI] Finally, we must choose values of Vao that encompass 
a wide range of possible physical values. 



4-1.1 Characterizing Finite Separations 

Strictly, our numerical simulations cannot treat infinite dis- 
tances. Therefore, we must propagate forward the 2-body 
star-star hyperbolic solution until the mutual distance be- 
tween the stars, r(t), reaches a specified distance r(t — 
tstart) = |7"start|. At this separation, the initial conditions 
for the numerical simulations are established. For any finite 
|?*start|; the components of the stars' velocity both parallel 
and perpendicular to the impact parameter segment will be 
different than their initial values at an infinite mutual sepa- 
ration. Denote the position and velocity components parallel 
to & by 2; and x, and the perpendicular components by y and 
y. Given values of n,Voo,b, and |r s tart| of the relative orbit, 
we seek to derive istart , 2/start , istart and y Bt ^t ■ 

Denote the hyperbolic anomaly as Eh- Then we adopt 
the same definition of r as in Pg. 45 of El <|l985h and Pg. 
85 of [Royl (|2005l ) such that 



cosh Eh 



1 + 



(24) 



Note that by this convention, r s tart is negative. Through ge- 
ometry, we have: 



= a h (e h - coshE h ) 



= V e| - lV cosh 2 E h - 1. 
Differentiating these gives 

dE h 



-ahV cosh 2 Eh — 1 



y s tart = a h \je\ — \ cosh E h 
yielding a total velocity of 

rf7" s tart 



dt 
dEh 
dt 



(25) 
(26) 



(27) 
(28) 



dt 



dEh 
dt 



±a h \/ e\ cosh 2 E h - 1 



7"start ah 



(29) 



where the RHS is due to the properties of a two-body hyper- 
bolic orbit. Thus, 



dE h 
dt 



ah ) 



(30) 



a-h\Je\ cosh E h — 1 



where the signs indicate the two possible velocity vectors 
along the orbit. Because we model the systems approaching 
one another, we adopt the upper sign. Subsequently, we can 
insert Eq. (30]) into Eqs. (27]) and (28} and use Eqs. O, (2}, 
(I24[) and p0|) to express the starting Cartesian elements in 
terms of given parameters: 
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% start — 



J/start = — bVoc 



-b 2 V£ + r start (r 8t artV2> - 2/i) 



(31) 
(32) 



% start — 



/' 



^"start 



' r star t 


+ 2/i 


^"start 


- 2/i 



-6 2 V£ + r start (r at a r tK2 - 2/x) 



2p 



J/start 



^start 

v/r Btart (r st artV^ - 2 M ) (6 2 V4"T7?) ' 



(33) 



(34) 



Finally, we convert these elements into the center of mass 
frame for the numerical simulation initial conditions. 

We wish to i) model an approximately equal approach 
and retreat for each simulation, and ii) sufficiently sample 
both the approach and retreat. Regarding i), in the reduced 
two-body hyperbolic problem, we need the time the approach- 
ing system takes to change \Eh\ to — \Eh\, or instead j smEh\ 
to — I sin Eh I : 



integrate 



= 2 T 



sinh E. 



= 2'- 



ItanhE 



dt 



= 2 



/ r start \/2fj, — rstart V"^ 



'"start Voo 



w 



VJ, + 2r a tartMl 



'start 



^ '"start 



(35) 



where r st art < 0. 

Regarding ii), suppose the longer of the planetary orbital 
periods is denoted by Ty. One wishes to sample a of these 
periods before the close encounter. Then, tintegrate/Tfc/ ^ 2a, 
or 



rstart 



-2W-^Vc< 



(36) 



We set r s tart from Eq. (|36[) . and set a — 1.2, in all of our nu- 
merical simulations in order to sample at least one planetary 
orbit both before and after the encounter. 



4-1-2 System Orientations 

As argued in Section 1, there is no apparent preferred di- 
rection for planetary system close encounters with respect to 
the Galactic Centre nor with one another. Therefore, we ran- 
domly orient both planets with respect to their parent stars 
and randomly orient both planetary systems with respect to 
one another. 



4-1.3 Impact Parameter Range 

Following previous work fe.g. lFregeau et al we express 

q as a multiple of the sum of the initial separations of both 



planet systems, so that q = p (aio + 020)1 where P is a con- 
stant. Using this form of the pericenter, we obtain from Eq. 
©: 

frmax " 



[P (aio + a.20)} 2 + [P (a w + a 2 o)], 



which is bounded from below as 
min(6 m ax)= lim b = P (a w + a 20 ) ■ 



(37) 



(38) 



The general expression for the upper bound is long, but may 
be simplified in specific cases. If we assume M s i = M S 2, 
M p i = M P 2 and aio = 120, which are the same assumptions 
adopted in our analytical cases, then 



b' = 2a20 



\ 



P 



P + 2 



Kir 



(39) 



where 6' denotes the value of 6 under the above assumptions. 
Hence, 



max (b' m 



lim b' 



= 2ai 



P[P + 



Ms2 + M p2 



2M, 



(40) 



which shows that the impact parameter may be arbitrarily 
large for a small enough planetary mass. 

For each set of simulations, we wish to sample a rep- 
resentative range of impact parameters. We chose P = 2.5 
and select values of b from out to o ma x according to b — 
VRAND x b^ax, where RAND is a low-discrepancy quasi- 
random Niederreiter number between zero and unity. The im- 
pact parameter is thus sampled according to its probability, 
and no weighting of the scattering experiment outcomes is 
necessary to account for the larger frequency of wide encoun- 
ters compared to nearly head-on encounters. 



4-1-4 6 an d V 

We choose three values of the planetary semimajor axis ratio 
(= e = aio/<i20 = 1, 10, 100), which represents a wide variety 
of both already observed planetary systems and systems with 
wide-orbit planets which have not yet been observed. In order 
to determine a range of plausible r/ values, reconsider Fig. 
[3] For aio = 1000 AU, our three values of e, and plausible 
velocity values in the field (10 km/s ^ Vco ^ 100 km/s), we 
obtain r] ranges of roughly [330 - 3300] , [120 - 1200], and 
[50 — 500], respectively. However, we need not restrict our 
r\ ranges to field values. We can also include typically slow 
cluster velocities of w 1 km/s (by reducing the lower bounds 
on the ranges by an order of magnitude) and hyper-velocity 
stars (by increasing the upper bounds by a factor of a few). 
Therefore, for each value of e, we choose between 15-20 values 
of 7] based on these broad ranges. 



4-1.5 Numerical Code 

We use a modified version of Piet Hut's fourth-order Hermite 
integratoi0, which we call SuperHermite. SuperHermite was 



That code is available at http://www.artcompsci.org 
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introduced in iMoeckel fc Verasl <|2012h . where the details of 
the implementation and code verification can be found. The 
SuperHermite code utilizes a P(EC) n method l|Kokubo et al.l 
1 19981 ) to achieve implicit time-symmetry. There is no pre- 
ferred dominant force or geometry (such as one central star 
or a circumbinary system). SuperHermite also contains col- 
lision detection; in all simulations, we set each star's radius 
to be a Solar radius and each planet's radius to be Jupiter's 
radius. Although SuperHermite can accurately determine ec- 
centricity variations many orders of magnitude smaller than 
observationally detectable values, in this study we consider 
only T > 10" 4 . 

As a safety measure, we choose a maximum allowable 
timestep for our integrations: 



^step,] 



2-7T7 



(41) 



Here, 7 is the fraction of the innermost orbit that the simula- 
tion is allowed to use as a timestep. We use 7 = 1/20. How- 
ever, SuperHermite's timestep choice will almost certainly be 
more conservative than this in all of our single-encounter sim- 
ulations. 



4-1.6 Other Considerations 

For each pair (e, rj), we ran iVtotal = 10 4 scattering experi- 
ments. Strictly, the imposition of a finite separation means 
that due to the center of mass frame shift, numerically 
e k o(t = istart) and a k o(t = tstart) deviate from their given 
initial values at infinite separations by a factor of roughly a 
few S. For example, eko{t = tstart) ~ 0.003. As we are pri- 
marily concerned with the change in eccentricity, this initial 
small nonzero eccentricity is not of concern. Regarding the 
change of orbital parameters, a*,(t > tstart) and et(t > tstart) 
do vary slightly as the stars approach each other, well before 
the close encounter. This variation, which increases with de- 
creasing |r a tart|, is natural and unavoidable, and is orders of 
magnitude less than the variation due to the close encounter. 



4.2 Simulation Results 



We present our cross-sections in Figs. [Tlllll Each successive 
figure shows a higher value of T, culminating with ejec- 
tion (T = 1). Each figure contains two panels; the left is 
for I Aei I and the right is for | Ae2 1 ■ The black circles, or- 
ange squares and blue diamonds respectively show the cases 
e = {1, 10, 100}. Each data point has vertical error bars; in 
some cases these are so small that they are not discernible. 

Due to symmetry, the black circles on both panels in each 
plot should be, and are, roughly equivalent. For most of the 
cross sections in the ejection figure (Fig. QTJ, j us t one data 
point was obtained for a particular (77, e) pair. Nevertheless, 
the plot demonstrates that ejection can occur, and predomi- 
nantly for the widest orbit planets. 

One perhaps surprising trend that is apparent in the left 
panels of Figs. [7][8] is that the normalized cross sections do 
not appear to be monotonic functions of e. Now we show 
how this trend indeed may arise naturally through analytic 
considerations. 



4-2.1 Semimajor Axis Dependence Explanation 

We reconsider the impulse approximation and the close con- 
figuration. Now remove the assumption aio = <J20- We seek 
the perturbation on Planet #1, whose initial semimajor axis 
is fixed, while the initial e is allowed to vary across simula- 
tions. This perturbation should be equal to 



|AVx| P i = \AV±\sis2 + \AV ± \ p2sl 

AVj_|s2pl + |AVj_|plp2 



(42) 



A similar analysis to that from Section [3] and the Appendix 
yields more complex formulae because here e / 1. We find 
that the resulting eccentricity excitation can be well approx- 
imated by: 



»(<0 
-if 



V2a w (2a w + e(b-2q)) 



bn (aio - eq) y/SjT+7) 



y/2a w [e (1 + e) bSr, 2 - 2a 10 (S V 2 + e (l + S + 8e 2 ))] 



brj^/5 (1 + e) [b8e V 2 (1 + e) - aio (5r? 2 + e (1 + S + 5e 2 ))} 



(43) 



We plot Eq. (|4*51) in Figs. [T2] and \TS\ in a regime which 
showcases the nonmonotinicity of efl (e). Note in particular 
how the orange curves are higher than the black curves, just 
as in the cross section plots. 

Further, the cross section itself is dependent on e: 



(aio + a 2 o) 
oc /? 2 + 2/3f%^ 



aio 
2 



2 , lN -2 
1 + - 



1 + 1 



-1 



(44) 
(45) 



where the constant of proportionality could itself be a com- 
plex function of e given, for example, Eq. (|43[) . 

We caution that these results are based on a single orbital 
configuration set and assume that the impulse approxima- 
tion holds, which becomes increasingly unlikely as the inner 
planet's semimajor axis is decreased (Eq. [7J. Nevertheless, 
they demonstrate how the dependence may be explained. 



4.3 Eccentricity Excitation Frequencies 

Having obtained cross sections, we can now determine the fre- 
quency with which planets' eccentricities are excited to par- 
ticular values. In particular, we are interested in the number 
of times over a main sequence lifetime that |Aei| > T occurs. 
Let us denote this number by A/", and the space density of a 
patch of the Milky Way as n and the main sequence lifetime 
as tMS (as in Section 1). Then 



a(|Aei| > T,e, »7)n 1 l4otMS 
ff(|Aei| > T,e, n)nr]tMS 
aio 



1000 AU 



0.5pc" 



2GM sl 5(l + e) 
aio V 1 + S 

tMS 

10 10 yr 



(46) 



where 
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Figure 7. Normalized cross sections for outcomes corresponding to |Aei| > 10 4 (left panel) and |Ae2| > 10 4 (right panel). 
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Figure 8. Normalized cross sections for outcomes corresponding to |Aei| > 10 3 (left panel) and |Ae2| > 10 3 (right panel). 



I = 0.016a norm (| Aei| > T, e, r,) r + J' * (47) 

is determined entirely from the numerical simulations, and is 
the only quantity in Eq. (|46[) determined by the numerical 
simulations. By setting J\f = f , one could obtain a "fiducial" 
value of A/" for wide orbit planets in the field with typical main 
sequence lifetimes. Because Af oc a% 2 , planets on tight orbits 
are well-protected. Nevertheless, a nonzero fraction of these 
planets will be affected. Thus, even if just a few percent of the 
~ 10 11 Milky Way stars host planets, millions of tight-orbit 
planets may be affected. 

Using our cross sections, we plot f , which represent fidu- 
cial values of N , in Figs. 1 1541 191 On the right axes of these 
plots, we indicate what the value of N would be for aio = 10 
AU, which is 10 -3 less than the aio = 1000 AU case. Further, 
we shade three regions on each plot, the horizontal extent of 
which correspond to lOkm/s ^ Voo ^ lOOkm/s for aio = 1000 
AU (left panels) and aio = 10 AU (right panels). The ver- 
tical extents have no physical meaning and were chosen to 
be nonintrusive. This velocity range corresponds to the typi- 
cal range of stellar velocities in the Galactic Disc. Therefore, 
these figures allow us to read off directly quantities of interest. 

Before analyzing the consequences of these plots, we first 



attempt to explain the trends observed in these figures based 
on our analytics. 

4-3.1 Explanation for the Frequency Trends 

The features in Figs. H5H19l are highly dependent on the range 
of b chosen. For example, if one chose b > (&stat >) (Eq- IA16]) 
exclusively, then the maximum eccentricity variation would be 
e oxt,max (Eg . I A24] ) . If instead b was sampled only inside &eje2> 
then the minimum eccentricity variation would be given by 
e$(f>) (Eq.EUl. 

Although the numerical integrations model interactions 
at random orientations, we can use our two limiting analytical 
cases in order to help explain the bumps in Figs. I15H18I In the 
far case, the eccentricity excitation monotonically decreases 
with both b and Voo (Fig. [5]). Alternatively, in the close case 
(Fig. [5]), qualitative differences are more varied. In that figure, 
we can superimpose horizontal lines which would be related 
to the values of T chosen in Figs. 1 1541181 Then we can count 
the number of instances when the Fig. [S] curves are above the 
horizontal lines: this yields a subset of Af '. 

Low horizontal lines, under e^ t min , all attain the same 

contribution for b < (6^t >)• However, for b > (b^t >)> l ne 
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contribution steadily increases as Voo is increased. This be- 
haviour is seen in the blue points of Fig. 1151 When the blue 
points tail off, effectively Voo has become high enough that 
the entire curve to the right of > i s under the horizontal 
line. The two rightmost blue points correspond to a value of 
Voo so high that the region b < b^ e 2 now becomes important. 
In this region, the dipping Fig. [5] curves dip low enough so 
that e^' t min < T, causing a large drop in M ■ This oscilla- 
tory behaviour is repeated in the blue curves of Figs. I16H18I 
as T, and hence the horizontal line in Fig. [S] steadily moves 
upward. By Fig. 1181 T is so high that increasing Voo serves 
only to decrease N ■ Figure ITTl shows the greatest detail, with 
a clear upward and downward trend plus modulations. These 
modulations likely naturally result from the random orienta- 
tions sampled, as it is important to recall that Fig. [6] models 
only a single, ideal configuration. 

Figure [19] is based on small number statistics and hence 
must be treated with caution. In particular, the linear upward 
trends of the black dots and red squares from the minima ex- 
plained in the last paragraph are based on single data point 
statistics. As evidenced by Fig. 1111 these single data points all 
have the same normalized cross section, meaning that M oar) 
(Eq. 146 p . The figure does demonstrate that ejection is pos- 



sible, even with stellar velocities towards the upper end of 
typical Disc velocities. 

Also, in all cases, despite the variations, N — > as Voo — > 
oo (Eqs. l46ll22l and l37|) . 



5 INTERPRETATION OF RESULTS 

Our results demonstrate that exoplanets are in constant dan- 
ger of losing their primordial eccentricities. The extent to 
which these eccentricities vary over time may i) provide a link 
to formation theories from the currently observed middle-aged 
systems, ii) identify the most dynamically excited regions of 
the Milky Way, and iii) impose a significant lower bound for 
the eccentricity variation of wide-orbit planets. 

5.1 Links to Formation Theories 

Classical core accretion typically forms planets within several 
tens of AU of the parent stars. These planets may begin their 
lives on nearly circular orbits, particularly if they are born 
in isolation. Without additional planets in the system, and 
provided that the formed planet is far away enough from its 
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CLOSE Case Impulse Approximation 




Figure 12. The change in e\ as a function of r\ for a few dif- 
ferent curves of constant e = a\Q /a^o in the close case. The plot 
demonstrates that orbital eccentricity variations are not necessar- 
ily monotonic functions of e. Consequently, neither are the cross 
sections, as shown in Figs. 17181 
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parent star to avoid tidal circularization, this planet can be 
perturbed only by external forces. 

Hence, the nonzero eccentricities of isolated planets may 
arise from planetary system flybys. The cumulative effect of 
planetary system flybys over a main sequence lifetime could 
eliminate the near-circular signature of a formation pathway. 
A few percent of planets on tight orbits could have their ec- 
centricities perturbed by 10 -3 (Fig. I16|l . which is comparable 
to the smallest observational e rrors yet achieved on plane- 
tary eccentricity measurements (|Wolszczan| [l994: Wels h et alj 
2012). A nonzero fraction of tight-orbit planets will experi- 
ence greater perturbations, with ~ 0.01% receiving a kick 
of over 0.1 (Fig. I18|) . Given that the current number of ob- 
served exoplanets is ~ 10 3 , we have not yet observed enough 
exoplanets, on average, to have detected a tight-orbit planet 
with such a large kick. Nevertheless, at least millions of such 
planets should exist in the Mil ky Way, given the re cent total 
exoplanet population estimate (|Cassan et alj|20l3 ). 

In multiple planet systems, the effect of flybys may be 
more pronounced. Two planets on the verge of dynamical 
instability could be driven to scatter off of one another after 
a sufficiently strong nudge from a flyby. Small, flyby-induced 
changes to particular dynamical signatures of formation, such 
as the circulation of the apsidal angle between two planets, 
can propagate over several Gyr to obscure the formative value. 

If, however, core-accreted planets are not born in isola- 
tion, and instead are continually perturbed in dense clusters, 
then these planets will attain a nonzero eccentricity. A de- 
tailed comparison of the relative contributions to a planet's 
dynamical history from its birth cluster versus its middle-aged 
interactions in the Galactic Disc may be crucial in determin- 
ing the types of planetary orbits seen in different Galactic 
environments. This study is a first step towards such a com- 
parison. Subsequent studies could focus on merging the two 
types of simulations, or at least consider fast interactions with 
initial non-zero planetary eccentricities. 

Nevertheless, we can provide a broad estimate here by 
com puting an eccentri city distribution similar to that of Fig. 
9 in iBolev et all (|2012r i, which illustrates a post-cluster plan- 
etary eccentricity distribution. Assume a population of any 
number of Solar-mass stars each with a 10 Gyr Main Se- 
quence lifetime and a space density of n = 0.5 pc -3 . Each 
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Figure 14. The minimum eccentricities achieved by the fraction of 
systems in 9 different model populations. All populations assume a 
Jupiter-mass planet is orbiting a Solar-mass star on a circular orbit 
at the same semimajor axis, a space density of n = 0.5 pc -3 and a 
main sequence lifetime of 10 Gyr. The symbols and lines correspond 
to: open black squares (Voo = 10 — 100 km/s, a = 1000 AU), filled 
vermilion squares (Voo = 10 — 100 km/s, a = 100 AU), filled right- 
pointing blue triangles (Voo = 10 — 100 km/s, a = 10 AU), open 
gray circles (V"oo = 10 — 30 km/s, a = 1000 AU), filled aqua circles 
(Voo = 10-30 km/s, a = 100 AU), filled left-pointing orange 
triangles (Voo = 10 — 30 km/s, a = 10 AU), open green diamonds 
(Voo = 80-100 km/s, a = 1000 AU), filled purple diamonds (Voo = 
80 — 100 km/s, a = 100 AU), and filled upward-pointing yellow 
triangles (Voo = 80 - 100 km/s, a = 10 AU). The plot broadly 
suggests regions of phase space where eccentricity excitation due 
to fast "middle-aged" Galactic Disc encounters may be comparable 
to or smaller than those achieved from other eras of the planet's 
lifetime. 



star has a Jupiter-mass planet orbiting on a circular orbit all 
of the same semimajor axis. Then we can use Eq. (|46[) with 
e = 1 and a given velocity distribution of the stars to compute 
M. Because <7 n orm(| Aei| > T, e, r/) is a discrete function com- 
puted numerically, we create an interpolating function based 
on those data points, for a given T. Then, for the velocity 
range of interest, we use the mean value theorem on the in- 
terpolating function to compute an averaged value of M for 
a given T. If M ^ 1, then we say that 100% of that popula- 
tion suffered an eccentricity of at least T. Recall that because 
our numerical integrations modeled single encounters, we do 
not know how additive the eccentricities are due to repeated 
perturbations when M 2. 

We apply this procedure to 9 populations: for a — 10, 100, 
and 1000 AU, and for flat distributions of Voo in the ranges 
[10-100 km/s], [10-30 km/s] and [80-100 k m/s]. The results are 
presented in Fig. 1141 A comparison with iBolev et all (2012) 
is difficult because of the different setups of the two papers. 
However, our Fig. [l4]does illustrate that a few percent of the 
planets at a — 100 AU experience eccentricity changes of at 
least 0.1, which may be comparable to those achieved in birth 
clusters. Further, for fast perturbers, the eccentricity change 
is a strong function of semimajor axis and a weak function of 
the velocity range chosen. Additionally, the values in the plot 
are dependent on n in a linear fashion, such that for dense 
environments, the fractions may increase by a factor of 2-3. 
Overall, more detailed comparisons are necessary. 



5.2 Dynamically Excited Galactic Regions 

Given that N oc nrja(rj) (Eq. 1461) . the extent of the plane- 
tary orbital disruption may significantly depend on the Galac- 
tic environment of the host star. Additionally, the migration 
history of an exoplanet host star through regions of differ- 
ing spatial density and velocity dispersions will affect the re- 
sulting perturbations on orbiting planets. Unlike planetary 
motion, stellar orbits are typically not closed, and can sud- 
denly transition from, for example, being ensconced in a dense 
tidal tail to traveling in the sparse region between two spiral 
arms. Even the region exterior to the Galactic Disc is com- 
plex: the two broadly overlapping structural components of 
the Milky Way's ha lo feature distinct spatial density profiles 
l|Carollo et alj|2007t ). 

Exoplanet host star velocities can vary by approximately 
two orders of magnitude. HIP 13044, which is thought to be 
of extragalactic origin, has a m easured systematic ve locity of 
300 km/ s with respe c t to S un (|Setiawan et al.ll2010t ). Alter- 
natively, iHelmi et al.l (|l999l ) suggest that relic debris streams 
from Milky Way formation show internal velocity dispersions 
of just a few km/s. In the Solar neighborhood (within 30 pc 
of the Sun), no single stellar vel ocity component exceeds 50 
km/s (|Nakaiima fc Morindl2012r i. 

These results suggest that the value of M may vary by 
a few orders of magnitude depending on the region studied. 
This variation may be important for characterizing the abun- 
dance and location of exoplanets when assessing the Milky 
Way's global population. In particular, for dense enough envi- 
ronments, wide orbit planets may survive only for a small frac- 
tion of the host's main sequence lifetime. Conversely, sparse 
environments would allow planetary systems to retain their 
formation signatures for several Gyr. Generally, regions closer 
to the Galactic Centre are denser, and hence perhaps har- 
bor more dynamically excited exoplanets than in the Solar 
neighborhood. Independently, this conclusion also arises from 
modelling the effect of galactic tides on exoplanets. 



5.3 Consequences for Wide-orbit Planets 

At least three exoplanets have been detected orbiting 
their parent stars at semimajor axes exceeding 10 3 AU 
llGoldman et all l201fj : iKuzuhara et al.l l201ll : iLuhman et~ail 
l201ll ). Several others are thought to orbit at separations 
of 10 2 AU - 10 3 AU. The population of wide orbit planets 
may be large, but remains difficult to distinguis h from the 
purportedly vast free-float ing planet population ( Su mi et al.l 
l201ll : iBennett et alj |2012| ) . Unlikely to have formed in their 
current locations via core accretion, wide-orbit planets per- 
haps alr eady represent the victims of inte rnal dynamical 
jostling (|Veras e t al. 2009; B olev et "ail 20121 ) or recaptured 
free-floaters ( Perets fc Kouwenhovenl 20121 ). Regardless, at 
these distances, these planets become even more susceptible 
to influence from external flybys. 

A wide orbit planet will typically have its eccentricity 
kicked by at least 0.1 roughly once over its host star's main 
sequence lifetime (Fig. I18[) . Further, the planet is likely to 
experience hundreds of kicks at the 10" 4 level (Fig. [IS). The 
probability of ejection is on the order of a few percent (Fig. 
I19[) . These values can vary by a factor of a few depending on 
the size of the intruder system's planetary orbit (e). Thus, 
wide orbit planets could represent an additional source of the 
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free-floating planet population , which cannot be explain ed by 
planet-planet scattering alone (|Veras fc Ravmondll2012l ') . Fur- 
ther, the significant eccentricity and semimajor axis kick given 
to wide orbit planets during their parent star's main sequence 
could has ten escape during that star's post-m ain sequence 
evolution jVeras et al.ll201ll ; IVeras fc Touj|2012l l. 



6 DISCUSSION 

Here we consider a few potential extensions to this work. 
First, we remove the assumption that the planetary masses 
are equal and estimate how our results might change. Second, 
we discuss other related few-body interactions in the Galactic 
Disc. 



6.1 Unequal Planetary Masses 

Here we briefly consider how orbital parameters might be per- 
turbed when the planetary masses are unequal. The ratio of 
planetary masses should be most important when the plan- 
ets are near each other during the close encounter of the two 
systems. Therefore, let us consider the close configuration, 
and specifically focus on the region around the planet-planet 
collision point (roughly bounded by 6 s tat,< and 6 s tat,>) 

Denote Si = M p \/M s \ and 82 = M P 2/M S 2- For simplicity, 
choose M s i = M S 2 and dio = 020, as in Section [3] If we 
carry out the same analytic procedure in that section and the 
Appendix, then we find that the eccentricity change of Planet 
#2 is similarly described by Eq. (JXTJ) , except now: 



2GAf.iV* (6aao-6(3-2* 3 )) 



(48) 



Zr, 



VL (4fl2o + b 2 (1 - 2<5 2 ) - ba 20 (4 + & - 3<5 2 )) (49) 



The result is the critical points around the planet-planet 
collision region become: 



"2(1 



('211 



4 - ygj + Sx (8 - 6S2) + 82 (8 + 9fe) 
2 - 4<5 2 

4 + y/5 2 + <Si (8 - 6S 2 ) + 82 (8 + 9S 2 ) 
2-AS2 



(50) 
(51) 



which is equivalent to Eqs. (|A15|) - (|A16|) when 81 = 82- 

We plot these critical points as functions of the mass 
ratios in Fig. 1201 The plot demonstrates that given a Jupiter- 
mass Planet #1, the region of planet-planet gravitational in- 
fluence changes by « 0.1a2o if Planet #2's mass is an Earth- 
mass versus a Jupiter-mass. In the latter case, the region of 
influence is greater. Note also how the asymmetry of the two 
critical 6 s tat points is enhanced when the planetary masses 
approach the stellar masses. 



6.2 Other System Configurations 

Scattering simulations for different hierarchical configurations 
of 4 bodies, or for more than 4 bodies, would provide a more 
complete picture of planetary orbital excitation from passing 
stars during the host star's middle age. However, the phase 
space to be explored is prohibitive. Nevertheless, because the 



few-body problem admits few analytical solutions, studies of- 
ten must rely on numerical integrations. 

Alternatively, in the Galactic Disc, the impulse formal- 
ism may be generalized to any number of bodies in any ori- 
entations. Although the resulting analytical formulas are un- 
likely to be as compact as those presented here, they - sub- 
ject to the assumption in Eq. (J7]l - would be able to sam- 
ple the entire phase space. Such a formalism could be use- 
ful, for example, in modeling how secular or resonant evo- 
lution of_jnulti^plaiwt_Sjrstem change naturally over 
time. IZakamska fc Tremaind (|2004l ) consider secular eccen- 
tricity propagation. For resonant systems, this same prop- 
agation might kick planets into a deeper or shallower mean 
motion resonance, if not out of the resonance entirely. Results 
from the Kepler mission i llustrate that there i s an abundance 
of ne ar-resonant planets (Li ssauer et al.ll201ll : iFabrvckv et alj 
120121 1, 

In cases other than the close and far cases, impulses 
would cause both a perpendicular and parallel kick. The net 
effect could be modelled as a single impulse. If a planet is 
on an eccentric orbit before the kick, then the true anomaly 
of the planet must be taken into account. In principle, one 
could remove the error bars associated with Poisson counting 
statistics in Figs. 171 i"T1 bv generating those figures analytically, 
and then generalizing the figures with a given distribution of 
eccentricities. In this way, one can also quantify the preference 
of a planet's eccentricity to increase versus decrease given an 
initial nonzero value. 

This formalism should also work for hierarchical sys- 
tems: those with stars, planets and moons. The Solar Sys- 
tem demonstrates that moons typically orbit planets within 
half of a Hill radius. Further, a planet's Hill radius is pro- 
portional to its semimajor axis. Therefore, wide orbit planets 
with moon^f] could feature a widely spaced moon orbit, one 
which extends to several percent of the planet's semimajor 
axis. At a distance of 1000 AU, such an orbit would be com- 
parable to the Neptune-Sun separation, and hence could be 
disrupted by passing stars. 

Also, giv en the possible va st population of free-floating 
giant planets (|Sumi et al.ll20lil ). passing giant planets might 
be more common than passing stars. Then, the resulting 
binary-single interactions with a passing fr ee-floater and a 
plane tary system could become important (|Varvoglis et al.l 
|2012| ). If a giant free-floating planet of mass M p were to pass 
by a system with a planet of mass M p orbiting a star of 
mass M 3 , then the critical velocity of this configuration is 
^{28 + l)/(2 + 8) w 71% of the critical velocity of the tradi- 
tional stellar flyby binary-single scattering configuration. This 
reduction in Vrit is not enough to claim that the system will 
be completely ionized; hence, this situation may be treated in 
a similar impulse situation as this work. In the perhaps more 
exotic situation of two pairs of free-floating planet binaries 
suffering a close encounter, Vrit would be reduced from the 
traditional four-star encounter by a factor of l/y6 ~ 32. This 
reduction is significant enough that ionization would be much 
more likely in that case for typical Galactic field velocities. 



5 Wide-orbit planets scattered out to their current locations could 
have retained moons, whether the moons were formed in the cir- 
cumplanctary disc or were captured satellites. 



© XXXX RAS, MNRAS 000,[Q-23 



200 
150 

^ 100 

2 70 

3 50 

s 

g 30 

20 
15 



Number of Instances of |Aei|>10" 



4-body Disruption of Primordial Planet Signatures 17 
Number of Instances of |Ae2|>10~ 4 









♦ • ■ 




♦ 
• 

♦ 


■ 

■ 


01O=<*2O 
«10=10fl20 ' 


■ 

♦ 




a 10 =100a 2 o , 


♦ B 







0.200 




ft 1 CO 




ft 1 ftft 


< 


0.070 


o 


0.050 


|| 

c 




a 


0.030 


- 

c 




i. 


0.020 




0.015 





0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 

Logio(Ko/V cri t) 



200 
100 

30 

10 

3.0 

1.0 

0.3 

0.1 




a 10 =a 2 o 
fll0=10«20 

a 10 =100a 2 o 
. 



0.200 
0.100 

0.030 

0.010 

0.003 

0.001 

3xl0" 4 

lxlO" 4 



0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 

Logio(K,/Vcri t ) 



Figure 15. The number of times over a typical main sequence lifetime (10 Gyr) in a Galactic region with space density of 0.5 pc -3 that 
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Figure 19. Same as Fig, 1151 but for ejections. 



7 CONCLUSION 

We have modeled the close encounter of two single-planet ex- 
osystems in the Galactic Disc, which mimicks a common oc- 
curence during middle-aged planetary evolution. We obtained 
analytical formulae and numerical cross sections which may 
be useful for future population studies of exoplanets in specific 
regions of the Milky Way. The resulting change in orbital pa- 
rameters for wide-orbit (a w 100— 1000 AU) planets is signifi- 
cant (with a typical Ae of several hundredths to over a tenth) 
and potentially measurable, suggesting that these planets are 
highly unlikely to retain a static orbit during main sequence 
evolution. Although tight-orbit planets (with a < 10 AU) 
are more resistant to orbital changes, millions in the Milky 
Way will be affected, and lose their primordial orbital signa- 
tures. The most dynamically excited Milky Way exoplanets 
are likely to reside in the densest Galactic regions. 
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APPENDIX A: ADDITIONAL ANALYTICS 

This Appendix expounds upon the analytical results in Sec- 
tion [3] The resulting formulae may be applied generally to a 
particular exosystem of study, and contribute to our analytic 
understanding of the general four-body problem. 

We can relate the effective impact parameters to the clos- 
est approach distances of all of the planet-star combinations 
through geometry. We obtain: 



b s i P 2 = (g±aio) 



Opls2 



(q ± a 2 o) 



bVl 



bVi 



bpi P 2 = (<7±aio±a 2 o) 



1 + 



bV 2 



1 + 



bV^ 



2 \ 2 



(Al) 
(A2) 
(A3) 



where the upper and lower signs correspond to the far and 
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close cases, respectively. In the close is the value 

of b which gives b p i P 2 = 0. Eq. Q yields: 



'v ,,u ^ \ / [aw + a.20) ( aio + o 2 o + \ 



(A4) 



When b < 6 m i n , then b p i P 2 < and the systems cross 
orbits. For a low enough b (when 6 s ip2 < or b p i S 2 < 0), the 
stars directly pass through the region in-between the other 
system's planetary orbit. 

Now we impose the analytic simplification described in 
S ect ion \3 . 2 . 1 1 and apply the fiducial values from Section [3.2.21 
to accompany the analytics. 

To gain physical insight into the following situations, con- 
sider how the perpendicular impulses from Eqs. |SJ - (fTTj) tend 
toward zero for both b — > and b — > oo. Therefore, each im- 
pulse is maximized for a particular finite value of b. These 
values are given by: 



&crit,sls2 
^crit,slp2 



Qcrit,plp2 



v oc 
frcrit,pls2 



(1 + 8) GMsi 



± O20 



(5 + 8) GM S1 
V* 



± O20 



( SGMsi 



± O20 



(2 + 8) GM S1 



V 2 

1 X 



± O20 



where the upper and lower signs denote the far and close 
cases, respectively. For planetary systems, 8 is small and hence 
the expressions for fe C rit,ai P 2 and fe C rit,pip2 may be shortened. 
For the fiducial values we adopted above, 6 cr it, s i S 2 ~ 1-97 AU, 
very close to the collision point of the stars. Further, in the 

1998.03 AU, 
In the close 



far case, b[ f 2 . - w 997.04 AU and b U) 



Ly crit,slp2 ' crit,plp2 

which are both a few AU away from oio = 020 



case, o crit slp2 



1002.96 AU and b 



.(c) 

crit,plp2 



2001.97 AU. In 



this case, note further that when 5 = 0, then b^} it plp2 = b n 



otherwise, b' 



(c) 



! crit,plp2 



is slightly higher: 



b (c) 

u crit,plp2 



1 + 28 



Q20 
^min 



U cir 



2 + (2 + 8) 



Vrirc.O 



For our fiducial case, 6 



(A5) 



0.0019 AU = 2.8 x 10° 



u crit,plp2 "mm 

km, which is smaller than the radius of the Sun but not of 
any of the Solar System planets. 

These critical values interact with one another to produce 
the interesting dynamics below. We first consider the far case. 

In the limit of b — 0, the stars will collide and impart 
a large perpendicular kick to the planets. The kick will be 
greatest at b = 6 C rit,sis2, which is about only 0.2% of 020 for 
our fiducial case. Nevertheless, reductions of Eqs. (|13fl and 
(|14|l show that the planet will be ejected (02/ — > 00, e%f ^ 1) 



for all b from to 6^, 



where 



(,(/) 



O20 

2V X 



[2Vcirc,0 " Voo + VVaa (Vac + 12V circ ,o)] (A6) 



or, at about 118 AU for our fiducial case. This result shows 
how a passing star at ~ 100 AU can rip a bound planet off 



of another star, even if the passing star is at opposition with 
the other star's planet. 

The eccentricity and semimajor axis perturbations ap- 
proximated by Eqs. (|15[) - (|16[ ) are independent of planetary 
mass because, in this case the planets are always "far" 

2020 = 2000 AU > b[H plp2 ) from each other0 Similarly, be- 

( f ) ( f) 

cause 6r-g + 020 > 0crit,pi P 2 1 we should expect the eccentricity 
and semimajor axis distributions to be smooth functions of 

This formalism allows us to estimate the contribution to 
Planet #2's eccentricity variation from the potential of Planet 
#1 alone: 



,(/) 



|AU± 



pls2 



I A VI 



plp2 



|AVl|p 2 
8b 2 



(A7) 



(A8) 



Equation (|A8|) shows that the relative contribution from the 
planet is an increasing function of b. The maximum contribu- 
tion is 



26<5(& + 02o) - (6 + 2o 2 



,(/) 



= e^\b^ 



1-28 



(A9) 



showing that Planet #1 completely dominates the evolution 
when 5 = 1/3. Similarly, we can quantify the change in Planet 
#2's semimajor axis from Planet #1 alone: 



,(/) 



(o 2 / - O20 ) 



Sir 



Planet #1 Only 



„(/)" 
'2/ 



(02/ 



(bZ 



a 20j Total 
2 / ry \2 



,(/)%(/) 2 
-It e X 



(ZiZty - (Sbz 2 y 



where 
Zi = 

z 2 = 
z 3 = 



b + 020 
020 

2(6 + 2o 20 ) 



U cir 



(6 + 2020) -2*6(6 + 020). 



(AIO) 



(All) 

(A12) 
(A13) 



Because it is expressed as differences of squares, Eq. (]A10[) 
readily reveals the conditions that will zero out the planetary 
contribution. 

Similar to the eccentricity, the relative contribution to 
02/ from the planet is an increasing function of b. In the limit 

b — s> 00, 



,(/) 



1 - 28 



(A14) 



6 Equation doe s not reduce to Eqs. (9) or (11) of 

iHeggie fc Rasid <199d) because the assumptions used to derive 
those latter two formulas are different: encounters are assumed 
to be slow, the encounter trajectory is assumed to be parabolic, 
and the impact parameter is assumed to be large relative to the 
planet-star semimajor axes. Further, the power-law dependence is 
reported in terms of pericenter distance. 
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FAR Case Impulse Approximation 
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<5 2 /(1-2<5) 2 



4000 



FAR Case Impulse Approximation 




4000 



Figure Al. Eccentricity and semimajor axis variation of Planet 
#2 due to Planet #1 alone in the far case. In the upper panel, Voo 
is varied; the resulting differences are negligible. In the lower panel, 

( f) ( f) 

S is varied, illustrating the sensitive dependence of and on 
the planet-star mass ratio. For distant encounters of the most mas- 
sive exoplanets and low mass host stars, the planetary contribution 
can represent several percent of the overall contribution. 



again showing that Planet #1 completely dominates the evo- 



lution when 8 = 1/3. Comparing 4'i„ an d e^ ax suggests 
that intruding planets have a greater capacity to alter other 
planets' eccentricities than their semimajor axes. 

Figure [ATI graphically illustrates Planet #l's contribu- 
tion to the evolution of Planet #2 in the far case. The 

plots demonstrate the contrastingly weak and strong depen- 
( f ) ( f) 

dencies of eir ' and wd' on Voo and 8, respectively. Also, 



,(/) 



> o 



(/) 



always . For the most mas sive-possible exoplanets 
(w UMj - 16Mj ; ISpiegel et alj[201ll ) and impact parameters 
of a few thousand AU, the planetary contribution may reach 
10% of the overall contribution. 

Now we perform a similar analysis for the close case. 
Here, where the position vectors from each star to their or- 
biting planet point towards each other, the resulting orbital 
parameter evolution is a more complicated function of b. Fig- 
ures [()] and |X2] may be helpful guides for the following discus- 
sion. 

In the limiting case of b = 6 m i n , the planets collide with 
each other, and the kicks on each other have no perpendicular 
component (as can be seen in Eg. 1111) . Further, the perpen- 



dicular kicks from both stars would cancel out. Therefore, in 
this limit - if the collision could be neglected - the planets' 
orbital parameters would remain nearly unchanged. However, 
even a slight nonzero distance between the planets would pro- 
duce strong perpendicular kicks, much stronger than the kicks 
from the stars, and cause the planets to escape. This kick is 
highest at b^ it plp2 , which differs from 6 m i n on the order of a 
giant planet radius (Eg. IA5|) . Therefore, around the vicinity 
of & m in, the planets either collide with each other or escape. 

As b deviates from & m i n , eventually the kick contributions 
from both Star #1 and Planet #1 will cancel out the back 
reaction from Star #1 on Star #2. The result is that planet's 
orbital elements would remain unchanged. This situation oc- 
curs at: 



do) 



0,20 



020 



2 - y/8 (4 + 8) 



1-28 



2 + y/S (4 + 8) 



1-28 



(A15) 
(A16) 



where the subscripts < and > indicate that fo sta t is less than 
or greater than 6 m in- For our fiducial case, fe^t < ~ 1942 AU 
and 4 c it,> ~ 2066 AU. 

For b < & s tat,<, the perturbation on Planet #2 becomes 
high as Star #1 approaches. In the vicinity of b w 020, where 
Planet #2 collides with Star #1, the planet is either ejected 
or destroyed for & c jc,2 b ^ 6 c j c ,i, where 

020 



(C) , 
cjc.2 

.<?>. 
eje,l 



2Uoc 
2a 2 o 



[2K irC ,0 + Voo + \/Uoo (Voo - 12Kirc,0)] , (A17) 

(A18) 



.(c) 
°eje,2' 



For our fiducial case, b[ 

J<0 



eje,l 



1074 AU and 6^ 2 w 926 AU. 



For b < b c ^ 2 , the perturbations on Planet #2 stay high as 



the two stars approach each other (bottom panel of Fig. [4j. 
The planet will escape for b < fe 



,(c) 
eje,3' 



where 



6 (c) 

°eje,3 



020 
2Voo 



|^2Virc,0 



+ V* 



x/Voo (Voo - 12V circ , ) (A19) 



In our fiducial case, 6' 



,0) 
eje,3 



137 AU. In between b 



(c) 
eje,2 



and 



°cjl 31 the eccentricity and semimajor axis variations are min- 
imized at 



2-V2 



I12U 



1 + 2 



V cir 



(A20) 



or, a value of w 587 AU for our fiducial case. The resulting 
minimum eccentricity and semimajor axes values are 



=(<=) 



(6 + 4^2) 



Voire, 
Voo 



Q.20 



1 - ( 6 + 4V2 



V cir 



(A21) 
(A22) 



(c) 



or e cxt.min = U.Sbb and a^ >min = 1.155a 20 . Eqs. (EH]) 
and (|A22|) imply that Planet #2 cannot remain bound at 
any b < b[f c l if Voo < 11.7V c irc,o ~ 10.8km/s. Therefore, 
Fig. [6] does not feature a black solid curve (which represents 
Voo = 10 km/s) for b < 1074 AU in either panel. However, 
when using this critical relation, one should remember that 
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Figure A2. Same as Fig. [6] except for different curves of S instead 
of Voo . Higher values of S have a marked effect on the region where 
both planets suffer a close encounter with each other, in which the 
contribution from the parent stars is negligible. Some critical points 
not marked in Fig. [6] are marked here. 



SOLID: a 2 changes; DASHED: e 2 changes 




b/AU 

Figure A3. Eccentricity and semimajor axis variation of Planet 
#2 due to Planet #1 alone in the close case. Unlike in Fig. [ATI 
Planet #1 may dominate the evolution over a large region of impact 
parameter phase space. 



the impulse approximation starts to break down as Voo de- 
creases according to Eq. 0. 



Now let us consider b > b'. 



(c) 



In this regime, the planet- 



planet interaction becomes negligible, and Planet #2's evo- 
lution is dominated by |AVj_| s i S 2 and |AVj_|sip2- These im- 
pulses partially, but not completely, cancel each other out, 
and admit the greatest net perturbation on Planet #2 at 



(,(c) 



(2 + V2 



(121, 



Vcirc,0 
Voo 



(A23) 



or ~ 3421 AU, which correspondingly results in maximum 
eccentricity and semimajor axes values of 



.(«=) 



(6 - 4V2) 



Voi 



«20 



.-4-\/2 



v cir 



Vo< 



(A24) 
(A25) 



,(<=> 



0.011 and a 



Cc) 



= 1.00012a2o. Equations 



(|A24fl and (|A25|I importantly let us explore if Planet #2 could 



ever be ejected when b > 010 + «20- Ejection is possible only 
if Voo < 0.34Vci rc ,o w 0.32km/s. 

Finally, to complete our exploration of the impact pa- 
rameter phase space, as b — > 00, the orbital changes asymp- 
totically approach zero. 

Planet #1 plays a much larger role in altering the orbit 
of Planet #2 in the close case instead of the far case. This 
can be seen by the dependence of 62/ on 5, even though in 
route to the derivation of e2/, we neglected a higher order 
term due to S. The fully general case (Eq. I14|) makes no as- 
sumptions whatsoever about S, meaning that the formula is 
just as applicable to four stars as it is to two stars and two 
planets. Hence, we use Eq. (|14|) in order to plot Fig. lA21 which 
illustrates the dependence on S. 

The plot provides an effective region of planet-planet in- 
fluence, perhaps interpreted as the region where planet-planet 
gravitational focusing is important. For S = 0.1, this region 
is nearly as large as < b < b^ t < . Note however, that the 

orbital parameter variations for 6 < 6^ t < are nearly com- 
pletely independent of S. For b > bstat <' greater values of 
S have an overall weaker effect, because of the locations at 
which the forces partially cancel out one another. 

Now we can estimate what fraction of Planet #2's orbital 
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changes are due to Planet #1 alone. We find that e^f is equal 
to 

8b 2 Vl{b-a 2Q ) 



[2GM sl +V£(b- a 20 )] [4a| + & 2 (1 - 25) - 2a 20 b (2 - 5)] 



(A26) 



which takes the same limit of 6/(1 — 25) for b — s> oo as (Eq. 
\M\i . Similarly, the limit of a£ ] as b -> oo is 5 2 /(l - 25) 2 . We 
plot these contributions in Fig. IA3I The left and right panels 
illustrate the dependencies on Kxj and 8 respectively. Like in 
the far case, here and ajf are greatly sensitive to S. 
Unlike in the far case, there is a region of impact parameter 
phase space where Planet #l's contribution dominates the 
evolution. For S > 10~ 2 , the width of this region can extend 
beyond 10 3 AU. Even for b > 10 4 AU, the contribution due 
to massive exoplanets may still be a few percent of the overall 
contribution. 
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